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Abstract 

The silo discharge process is studied by molecular dynamics simulations. The development of 
the velocity profile and the probability density function for the displacements in the horizontal 
and vertical axis are obtained. The PDFs obtained at the beginning of the discharge reveal non- 
Gaussian statistics and super diffusive behaviors. When the stationary flow is developed, the PDFs 
at shorter temporal scales are non-Gaussian too. For big orifices a well defined transition between 
ballistic and diffusive regime is observed. In the case of a small outlet orifice, no well defined 
transition is observed. We use a nonlinear diffusion equation introduced in the framework of non- 
extensive thermodynamics in order to describe the movements of the grains. The solution of this 
equation gives a well defined relationship (7 = 2/(3— q)) between the anomalous diffusion exponent 
7 and the entropic parameter q introduced by the non-extensive formalism to fit the PDF of the 
fluctuations. 
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I. INTRODUCTION 



The process of gravity-driven silo discharge could naively be expected to be a simple 
one, yet we lack a well defined theoretical framework to explain the experimentally observed 
grain dynamics in terms of fundamental interactions. During the silo drainage, particles 
evolve under the action of gravity and interact between them through inelastic collisions 
in a complicated way, where intensive variables -such as density or temperature- are not 
well defined. Under these circumstances, the displacement and velocity fluctuations could 
be of the same order than its mean values, and the net force acting on each particle could 
even present strong deviations compared to the gravity force, therefore introducing unwanted 
effects like segregation or jamming . Two models were generally considered successful 

in explaining the global characteristics of the flow (the velocity profile) inside the silo. One 
of them is based on a continuous approach and the other takes into account the discrete 
nature of the media. The continuum model 4| uses concepts like elasto-plastic potentials 
introduced in the framework of continuum mechanics, in which the mean velocity field can 
be obtained. Nevertheless, this approach looses validity near the outlet orifice and does not 
bring any information about possible microscopic effects like mixing or segregation. The 
second one is based on the idea that each individual particle executes haphazard movements 
which can be treated as a random walk. This notion was originally introduced by Litwiniszyn 
in the sixties and is commonly called diffusive void model j^]. This model focuses on the 
movement of voids injected at the outlet of the silo. Such voids diffuse upwards, exchanging 
their positions with the grains, which in turn move towards the orifice at the bottom. In 
spite of their different premises, both approaches give essentially the same results for the 
mean velocity profile, but in the latter case the profile depends only on a single parameter, 
that can be regarded as a characteristic length a related to diffusion. 

In this work we present numerical simulations of a silo discharge process, including the 
beginning of the operation. We use three outlet diameters to study the behavior of the 
velocity profile, demonstrating an evolution between the transitory and stationary states. 
Results for the PDFs of the displacements of individual grains reveal non-Gaussian statistics 
and super-diffusive behavior at the beginning of the discharge. In agreement with exper- 
iments, non-Gaussian to Gaussian PDF's transition is observed in the stationary regime. 
Finally, we show that the complete sequence of dynamical states displayed by the particles 
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at the beginning of the discharge can be interpreted in the non-extensive statistical mechan- 
ics framework introduced by Tsallis {(|, which can be used to fit the obtained PDFs if the 
anomalous scaling of mean-square displacement as a function of time is taken into account. 

II. THE DIFFUSIVE APPROACH 

The diffusive void model yields a well defined relationship between the fluctuations of the 

n q 

positions of the particles and the characteristic length used to fit the velocity profile [5|, |7j . 
3ut recent experiments have evidenced some discrepancies concerning these predictions. In 
8|, a high-resolution particle-tracking experiment of a silo discharge is reported. The es- 
timated value for the parameter a -around 2 or 3 particle diameters- is much larger than 
the one predicted by the aforementioned authors. Choi et al. also reported the observation 
of non-Gaussian statistics and a super-diffusive regime of the grain displacements at suffi- 
ciently short temporal scales. (This temporal scale is defined by the time it takes for a bead 
to fall a distance of its own diameter). The probability density functions (PDFs) obtained 
for the particle displacement in the vertical and horizontal directions at this scale are both 
fat-tailed. For larger displacements both PDFs evolve toward a Gaussian shape. The grains 
undergo a transition from a super-diffusive to a diffusive regime in a coarse grained scale. At 
these scales, the velocity profile is stable and the diffusive models seem to regain their valid- 
ity. More recently, a similar experiment was reported where the non-Gaussian fluctuations 
remain even for the coarse grained scale 9]. 

Let us to go through the fundamental assumptions of the diffusive models. The line of 
reasoning proposed by Litwiniszyn later developed by Mullins basically introduces 
two fundamental assumptions: (i) the particle moves following a random walk in which the 
interaction between horizontal an vertical displacements can be neglected, and (ii) the time 
elapsed in one jump is almost the same than the one corresponding to a free fall. Assuming 
these hypothesis, in a typical experimental situation the diffusivity of the particle scales 
with its diameter d as D ~ \fgd A , and the parameter a -which Mullis named "jumping 
distance"- can be expressed as D/V = 3/8d, where V is the mean velocity in the vertical 
direction. The parameter a is relevant to check the theory against the experimental results 
because it is all that is needed to fit the velocity profile. The same parameter was derived 
from probabilistic arguments for an hexagonal array of particles js|, yielding the value dj 
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FIG. 1: Average vertical velocity profile inside the silo for orifices of size 11 (top) and 3.8 d (bottom) 
at different stages of its evolution. Time increases from left to right. The last picture corresponds 
to the fully developed regime. The localized structures correspond to a well defined region where 
strong velocity gradient exist. Note than this zones does not totally disappear in the stationary 
regime of the small outlet orifice. 

Taking into account the strongly biased character of the diffusive movements in the vertical 
direction, it is possible to write an analytic expression for the velocity profile (assuming a 
point source of voids) [3]: 



v (r, z) cx exp 



(1) 



x 2 + y 2 . 



where r 2 

The same profile was derived from phenomenological arguments by Nedderman & Tiiziin 
at the end of the the seventies |lo| . Following a kinematic reasoning, this model reproduces 
quite well the velocity profile inside the silo by assuming that the horizontal component of 



the velocity u is straightforwardly proportional to the gradient of the vertical velocity v: 



Considering the granular media as an incompressible material, an equation can be de- 
duced for the velocity profile which is isomorphic to the one obtained following the diffusive 
model: 

oz OX A 

In this kinematic model, the temporal variable has been replaced by the coordinate z. 
The solution to this parabolic equation corresponds to a Gaussian profile whose shape only 
depends on the parameter B: 

v(x,z) = ^2=e-^ B * (4) 
where Q stands for the volumetric flow. 

The parameter B appearing in Eq. [3] plays the same role than the "jumping length" a 
introduced by diffusive models. But contrary to them, the kinematic model does not provide 
any explanation about the dependence of B on the typical variables of the problem -such 
as particle diameter, silo size or any other characteristic length. In all the cases reported, 
the value of B lies between 2 or 3 times the particle diameter. Besides, B depends on the 
position inside the silo jll ]. 

Although the experimentally determined value of B does not match the jumping length 
predicted by the diffusive models either, both approaches rely on just a single parameter to 
describe the flow. This makes these kind of models appealing when trying to describe the 
movements of the particles in the discharge process. 

Recently, a more general description was introduced 1^] in order to provide an expla- 
nation for the discrepancy between the experimentally determined figures for the "diffusive 
length" and the values predicted by the diffusive approach. This new model assumes that 
the voids injected at the outlet of the silo do not match the size of just one particle. Instead, 
the injected void spreads through several grains, which will then move cooperatively. The 
resulting "cluster" , which spreads the void throughout a group of particles, is called a "spot" . 
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Spots move upward inside the silo and the grains affected by the spot carry out small move- 
ments toward the base. As the grains affected by the spot must move in a concerted fashion, 
it is not surprising that the displacements of nearby particles present strong correlations. 
The resulting spot diffusivity allows to recover the values for the jumping length obtained 
in experiments. 

Let us note that in all the cited references the discharge process must have reached a 
steady state regarding the flow, far enough from the beginning of the discharge. Moreover, 
the size of the orifice was large enough to avoid jamming events like those described in 
{lsI ]. Remarkably, the universality of the results may be lost when the size of the outlet 
orifice is small (see the case of small flow in reference 8j). Regardless of the inherent 
complexity of these states, a diffusive approach might be suitable to describe the particle 
displacements, although the simple assumptions introduced by the model could lose their 
validity. Early stages of the discharge process involve a densely packed media subjected 
to a strong, inhomogeneous shearing. Furthermore, even in the case of homogenous shear, 



densely packed granular media is prone to complex dynamics, as reported in 14 



where concepts like particle mobility and diffusivity to are invoked to describe the particle 
dynamics. 



III. NUMERICAL SIMULATIONS 



We have carried out molecular dynamics simulations of disks in two dimensions [17j. A 
simulation begins with 5000 disks arranged in a regular lattice; they are then given random 
velocities, which have a Gaussian distribution. The disks are allowed to fall under gravity 
through a conical silo. Below this hopper lies a flat bottomed silo where the grains are 
deposited. This is the preparation phase, which is purposely long and complex in order to 
break the correlations that the initial regular arrangement of the grains might induce in the 
dynamics. Once all the grains have fallen into the flat silo we wait until most of the kinetic 
energy is dissipated. Finally, we open an outlet at the bottom allowing the grains to fall 
and we start our measurements. 

The walls of the two silos used are built with grains. The interaction between grains is the 
same than the interaction between grains and walls, but the latter are fixed in their places. 
Thus the walls are rough and cause dissipative collisions. The width of the silo base is 50 
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grain diameters and the height reached by the grains when the silo is full is approximately 
twice that length. These dimensions guarantee that the flow rate does not depend on wall 
or filling effects. In addition, the filling height changes only slightly during a simulation. 

The model for the forces describing the interaction of two particles i, j consists of normal 
and transversal contacts as well as dissipative terms: 



where 



F n = - InVl, (5) 

F t = - min (fi\F n \, 7* l) ■ sign (y'Jj (6) 



v tj = • s + (ui + Uj) (7) 



Eq. (jSJ) gives the force in the normal direction of the contact. The first term is a restoring 
force proportional to the superposition £ of the disks. The 3/2 exponent arises from the Hertz 
theory of the contact. The second term is a dissipation proportional to the relative normal 
velocity of the interacting disks with damping coefficient j n . Eq. (j5J) is the sliding component 
of the damping force. It implements the Coulomb criterion with friction coefficient fi. The 
damping in the transverse direction is proportional to the shear velocity given by Eq. (JTj) 
and a transverse damping coefficient j s . In this equation s is a unit vector tangential to 
the disks at the contact point, and Ui and Uj are their angular velocities. Thus, in this 
scheme of forces we have a restoring force which prevents grains to interpenetrate (although 
a very small penetration is needed) along with damping terms in the normal and tangential 
directions which dissipate energy during the contact. This dissipation is among the most 
prominent characteristics of granular media. 



The values of the coefficients are, in reduced units, k n = 10 5 mg/d, 7 n = 100 mJg/d, 



7 S = 300mJg/d, and \i = 0.5. The integration time-step used is 1.25 _4 r with r = yd/g, 
and m, d and g stand respectively for the mass and diameter of the disks and the acceleration 
of gravity. 

The equations of motion were integrated using the velocity- Verlet algorithm and we used 



a neighbor list [18| to reduce the computational effort. 
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FIG. 2: Normalized PDFs for (a) horizontal and (b) vertical displacements. These distributions 
have been obtained in a temporal window where the averaged particle displacement is lower than 
a particle diameter. The dotted line is a Gaussian. 

A. Velocity profiles and probability density functions 

Using three different diameters of the exit orifice (namely, 3.8 d, 11 d and 16 d) we studied 
the evolution of the vertical velocity profile. The smaller and larger diameters where chosen 
because they belong clearly in two distinct regimes: for the former, the flow can be inter- 
mittent, while for the latter it is not l^ . We compute the velocity profile from the moment 
when the outlet is opened until the moment when the velocity profile becomes stationary. 
For each orifice size we perform 20 independent simulations and average them to obtain the 
final result. 



8 



1 

ft 





D-3J14 1 


□ 


D- lid 1 








6 10 



FIG. 3: Normalized PDFs for (a) horizontal and (b) vertical displacements. These distributions 
have been obtained in a temporal window where the averaged particle displacement is larger than 
a particle diameter. The dotted line is a Gaussian. 

In fig. [T] we show the evolution of the averaged velocity profile (only its vertical component) 
for the 3.8 and 11 d orifices. It is evident that in both cases groups of grains move downward 
together at the very beginning of the discharge, while structures that can be described as 
"bubbles" can be seen moving upward. Let us stress that these spatial structures do not 
correspond exactly to the spots introduced in [12] . These bubbles are zones where the mean 
velocity is larger than the bulk and its evolution reveals the intermittent regime at the 
beginning of the discharge. When these bubbles disappear, the characteristic stable flow 
profile is developed and the asymptotic velocity converges to a Gaussian profile. In the 
asymptotic regime (for the bigger orifices) the velocity profile in the middle of the silo can 
be fitted using a numerical "diffusivity length" of 2.2 ± 0.2 d which is in excellent agreement 
with experiments j§J, lit 
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FIG. 4: Mean squared displacements for the horizontal and vertical displacements in logarithmic 
scale. Circles correspond to the x component and triangles to the z component, (a) 16 d silo. The 
ballistic character for small displacements is obvious; once the displacement exceeds a distance of 
about the diameter of the particle the movement becomes diffusive, (b) 3.8 d silo. For small exit 
orifices it is difficult to distinguish the two behaviors and the crossover at about one d mentioned 
for the large orifice. The movement of the particles seems to be super diffusive at all scales. 

The results obtained with the smaller exit orifice are similar but an important difference 
arises: the characteristic time needed to reach a stable flow grows dramatically as the outlet 
diameter decreases. This feature is related to the fact that an increasing number of bubbles 
appear in the system, inducing an intermittent flow in the silo. Under these circumstances, 
strong velocity gradients can be observed near the exit orifice. This implies that some 
assembly of particles can have a relatively low velocity that may enable them to form an 
arch near the outlet, and stable jamming events can arise if the particles arrest the flow. 
This intermittent regime will be studied anywhere. 

We measured the displacements of individual grains during the discharge. In order to 
compare our results with those obtained experimentally, we chose a sampling time such that 
the mean displacement of the particles during it is of about 0.01 d (which is approximately 
the same than in the experiments). The particles were tracked in a window with different 
dimensions (10x10, 15x15 and 20x20 d) at the center of the silo. We obtain the same results 
regardless of window size. When computing the PDFs of the vertical displacements we 
subtract the corresponding component of the mean flow velocity. 

In fig. [2] we show the PDFs corresponding to the fully developed flow in semilogarithmic 
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scale. The displacements in each direction are normalized by their standard deviation. We 
see that the PDFs are essentially the same for the vertical and horizontal displacements. 
Besides, they are clearly non-Gaussian, with apparent differences both in the central region 
and in the tails, which are fat. As reported in |8(|, the fluctuations in the horizontal direction 
evolve toward a Gaussian profile when the displacements are larger than the diameter of 
the particle. On the other hand, the PDF corresponding to the fluctuations in the vertical 
direction remains long-tailed even for longer periods of time (see fig. 0b) wich is consistent 
with the experimental result reported by Moka et. at. {9]. 

In fig.HJwe show the mean squared displacements in each direction as a function of particle 
displacement (normalized by the particle diameter d) for a small and a large orifice. In the 
case of the largest orifice (fig. ffla), a well defined crossover between two different regimes 
is displayed. For displacements up to about one particle diameter, the variance displays a 
ballistic behavior, as in molecular fluid transport 19) . A normal diffusive regime can be 
used to describe particle fluctuations equal or larger than one diameter for big orifices. 

On the contrary, for the small orifice (fig. Hlb) there is not a well defined transition from 
ballistic to diffusive regime. Furthermore, subparticle displacements are superdiffusive. This 
is consistent with the fact that the corresponding PDF (fig. [3]) remains fat tailed at all scales. 
Remarkably, a subparticle superdiffusive regime was reported in the experiment 8| for all 
the outlet orifices studied. Such an effect would be associated to the strong influence of the 
lateral wall on the particle dynamics where the particle tracking is performed. 

It is remarkable that a single parameter like the one introduced by the kinematic model 

n Q 

[ICj or the jumping length introduced by Mullins [7j can be used to obtain the shape of 
the velocity profile, at least for large enough orifices. Unfortunately, the diffusive length 
B introduced by simple phenomenological arguments can not be easily associated with the 
microscopic grain dynamics. The diffusive model also predicts a Gaussian velocity profile 
that depends on a single parameter. Seemingly the meaning of this parameter is the same: 
a "characteristic diffusive length"; but care is needed in order to compare both two. The 
parameter introduced by Mullins is the ratio between diffusivity and the characteristic ve- 
locity of the grains (or the voids, for that case) during a "jump" comparable to its own 
diameter, whereas B is just a phenomenological parameter introduced to link the two ve- 
locity components (as no prediction is provided by the model, it must be obtained from 
experiments). But there is not any a priori a reason to consider both parameters equivalent 
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from a theoretical point of view. 

Admittedly, the predictions of the diffusive models (including the void model, which 
assumes the unrealistic situation of the voids executing a simple Bernoulli random walk in 
a regular lattice) do not match quantitatively the fitting parameter for the velocity profile 
found in numerical simulations or experiments. As noted above, a new model was introduced 
to tackle this discrepancy jl^ . This model proposes a cooperative mechanism or "cage 
effect" to represent the particle void interaction. This effect gives a correct estimation of 
the velocity fitting parameter but does not explain the microscopic diffusive regime. 

The spot model infringes one of the most important hypothesis of the diffusive models: 
void and particle displacement are no longer symmetrical. The symmetry plays a crucial 
role in the calculation of particle diffusivity by determining the time that a grain needs to 
migrate a distance equal to its own diameter. When a void enters the silo, the particle must 
perform a ballistic flight in a time scale comparable to the free fall time. The numerical 
simulation confirms the ballistic flight in the fully developed flow (Fig H]) but in a temporal 
scale considerable larger than a free fall. In our simulations the flight time is three or 
four times longer than the free fall time, depending on the position in the silo. This spatial 



llj that reveal a variation of the parameter 



dependency also agrees with experimental results 
B with the particular place of the silo where it is measured. Working is in progress in this 
line and results will be discussed elsewhere. 



IV. A NON-LINEAR DIFFUSIVE APPROACH FOR THE BEGINNING OF THE 
DISCHARGE 

Let us now describe the particle dynamics at the beginning of the discharge. This is 
important in order to understand the origin of jamming and the existence of a critical 
radius [3] beyond which the flow is never interrupted. In fig [5] the normalized PDFs for the 
displacements in each direction are displayed. They are clearly non-Gaussian and differ a 
little bit from the ones obtained for a stationary flow. This is indeed as expected. At the 
early stages of the discharge, the particle motions should be more correlated than afterwards. 
The material needs to dilate in order to flow, and the dilation takes place preferably along 
the regions where more empty volume is available. Under these circumstances, contact 
network or force chains can persist for a finite time and eventually block the outlet orifice. 
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FIG. 5: Normalized PDFs for (a) horizontal and (b) vertical displacements at the early stages of 
the discharge. These distributions have been obtained in a temporal window spanning from the 
beginning of the discharge to the moment when the mean displacement is twice the diameter of 
the beads. The dotted line is a Gaussian. The continuous line is Eq. [10] 



At the beginning of the discharge, the packing fraction will be lower than the corresponding 
to the ideal close packing of an hexagonal lattice and an important number of arches 20] 
will introduce long range interactions between particles. At the same time, each particle 
will be affected by non trivial stresses along the vertical and horizontal directions. 

The mean squared displacements scale approximately as (Ax 2 ) = (A2 2 ) oc t 4//3 for all the 
orifices studied (Fig. [6]). Moreover, there is not a clear crossover between superdiffusive to 
diffusive regime even for big orifices. These results agree with the experimental observation 
that even for large orifices is necessary to wait an easily mensurable time before reaching 
the stationary regime. Clearly this time strongly depends on the orifice size. 

The anomalous scaling and the PDFs obtained at the beginning of the discharge seem to 
be universal and suggest that it is necessary to introduce a generalized expression to describe 
the "walk" or displacements of the particles. 

Some theoretical models have been introduced to describe these familiar features of granu- 
lar media using alternative non standard thermodynamical descriptions, which resort to con- 
cepts like granular temperature 2lj. But these dynamical states are in fact quasi-stationary 
0, Q whereas in our case it would not be pertinent to apply a concept like temperature at 



the beginning of the discharge. Besides, recent numerical 16j and experimental works 
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FIG. 6: Mean squared displacements for the horizontal (a) and vertical (b) displacements at early 
stages of the discharge. The horizontal component shows clearly a superdiffusive behavior in all 
the cases. The vertical component displays a similar behavior although it is not as obvious. The 
curves plotted would become smoother if the number of averaged discharges is increased. 

have shown that when a system of grains is affected by shearing effects, anomalous scaling 
relationships are general. As argued by many authors, any transport mechanism where the 
variance of particle displacements scales with time as t 1 with 7 7^ 1 is connected to anoma- 



lous diffusion 



23l . |24| . A great variety of physical systems can display fluctuating transport 



mechanisms where such anomalous scaling can be found 251 ] . The paradigmatic example is 



23(. In such a case, randomness in the 



a random walk performed in a disordered network 
lattice induces inhomogeneous transition rates and as consequence the diffusive process is 
anomalous. In our simulations, the PDFs for the displacements indicate that the particle 
diffusion at the beginning of the discharge is anomalous. We will now introduce a formal 
framework to describe the anomalous behaviors. 

There are two important attempts to generalize the normal diffusive equation to include 
anomalous dynamics. The first one remains linear but introduces the concept of fractional 



derivative in its formulation 



261 ] . This formalism will be not discussed here due to the fact 



that it does not display a transition to norma! 



to switch to a nonlinear equation (NLDE) 27] ; in the simplest version it reads: 



diffusion at longer times. The second one is 



dp(x,t) _ D d 2 \p(x,t)\ 



(8) 



dt dx 2 

where v is a real number. One of the consequences of this functional dependence is that 
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the diffusivity has now a complex dependence on the density probability function p(x, f . 



The microscopic origin of such an equation has been discussed by many authors 



27 



28 



29|, 



specially in the framework of the Tsallis entropy formalism [2j, |25j. Such formalism was 
proposed to analyze long range interacting systems and its non-extensive properties. Since 
its introduction, Tsallis statistics has been used in a broad range of physical problems where 



non-Gaussian PDFs and anomalous diffusive behaviors appeared 
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framework of a granular gas 



33]. 



321 ]. even in the 



Tsallis formalism is important here because it gives an analytic solution to Eq. (JSJ) in- 
troducing a new parameter q called entropic parameter |24j. Considering a ID, case the 
solution for this equation is given by: 



p(x,t) 



-x 2 /[(3-g)(£>t) 2/(3 - 9) ] 



(9) 



v ^^( J Dt)i/(3_ Q ) i 

where q = 2 — u < 3. It has been shown [24| that the scaling between time and mean 
square distance can be expressed as a function of q through the expression (x 2 ) oc t 2 ^ 3 ^ q \ 
provided that the second moments of the displacement distribution remain finite, which is 
our case. 

Let us consider the case where (£2) — (z 2 ) oc t 4//3 (Fig. [6]). Comparing this numerically 
obtained relationship and the former equation, we can estimate the value for the q parameter: 
q = 3/2. With this value, we can write the analytical expression for the solution of Eq. ([8]) 
normalized by its standard deviation [25] as: 



p [X) 



2/tt 



:i + x 2 



(10) 



This function corresponds to the solid line plotted in fig. [21 The agreement between the 
numerical result and the analytical prediction indicates that a NLDE like Eq. [8] represents 
fairly well the evolution of the particles' fluctuations both in the vertical an horizontal 
directions at the early stages of the discharge. 



Furthermore, as time increases the so 
tend to a Gaussian for values of q < 5/3 



utions provided by the Tsallis formalism to Eq. 



251 ]. as observed in our simulations when the system 



reaches the stationary regime (for the x component). Although the meaning of the diffusivity 
constant in normal diffusion can be related to well defined concepts like temperature and 
mobility, the diffusivity introduced in Eq. (jSJ) lacks any obvious physical meaning. At this 
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stage of the discharge, the diffusivity D would be determined by the complex interactions 
between flowing particles. They give rise to a transport coefficient with a physical meaning 
unavoidably more complex than just a characteristic length obtained from a well defined 
particle velocity 



V. DISCUSSION 



We have studied the temporal fluctuations of particle displacements in the discharge of a 
silo by gravity. We show how at the early stages of the discharge the particle displacement 
presents distinctive features, as previously suggested by many authors. For that regime, it is 
possible to derive the PDFs for the displacement fluctuations of the particles from the non- 
extensive entropy theory introduced by Tsallis. We have shown that this formalism allows 
the treatment of diffusive systems presenting anomalous behavior through the introduction 
of a nonlinear diffusive equation. This equation would then provide a formal framework to 
represent the dynamical evolution of the PDFs from the beginning of the discharge to the 
fully developed regime. The introduction of this formalism might also offer an explanation 
for the jamming probability introduced in previous works [2I. [l3]. 

For the stationary regime, we corroborate the non-Gaussian features of the particle fluctu- 
ations for small displacements. Nevertheless, we found that the subparticle diffusive motions 
are ballistic, as is usual for the molecular fluids. In our opinion, the superdiffusive regime 
reported in some events Q cou.d b e r e,te d to th e effects of UrtenU w* on the p artlcl e 
mobility. It is a widely accepted fact that the typical boundary layer in granular material 
spans around 10 particle diameters, so the lateral size of the silo could cause strong effects 
on the particle displacement. 

The transition between ballistic to normal diffusive regime validates the applicability of 
diffusive models. The discrepancy observed between the characteristic length used to fit 
the velocity profile and the one predicted by the model is caused by a wrong estimation 
on the scale for the typical subparticle displacement velocity. The non-Gaussian nature for 
the PDFs for the vertical coordinate should be certainly linked with this scales and will be 
studied in the future. 



16 



VI. ACKNOWLEDGMENTS 



This work has been supported by project FIS2005-03881 (MEC, Spain), and PIUNA 
(University of Navarra). R.A. thanks Friends of the University de Navarra for a scholarship. 



[1] A. Samadani, A. Pradhan, & A. Kudrolli, Phys. Rev E 60, 7203 (1999). 
[2] K. To, P-Y Lai, & H. K. Pak, Phys. Rev. Lett. 86, 71 (2001). 

[3] I. Zuriguel, L. A. Pugnaloni, A. Garcimartm & D. Maza, Phys. Rev. E 68, 030301(R) (2003). 
[4] R. M. Nedderman, Statics and Kinematics of granular Materials, Cambridge Univ. Press, 

Cambridge UK, 1992. 
[5] J. Litwiniszyn, Bull Acad. Pol. Sci. 11, 593 (1963). 
[6] C. Tsallis, J. Stat. Phys. 52, 4579 (1998). 
[7] W.W. Mullins, J. Appl. Phys.,43, 665 (1972). 

[8] J. Choi, A. Kudrolli, R. R. Rosales, M. Z. Bazant, Phys. Rev. Lett. 92, 174301 (2004). 

[9] S. Moka, P. R. Nott, Phys. Rev. Lett. 95, 068003 (2005). 
[10] R. M. Nedderman, U. Tiiziin, Powder Technology 22, 243 (1979). 
[11] J. Choi, A. Kudrolli M. Z. Bazant, J. Phys-Cond. Matt. 17, S2533 (2005). 
[12] M.Z. Bazant, Mechanics of Materials 38, 717 (2006). 

[13] I. Zuriguel, A. Garcimartm, D. Maza, L. A. Pugnaloni & J. M. Pastor, Phys. Rev. E 71, 

051303 (2005). 
[14] H. Makse & J. Kurchan, Nature 145, 614 (2002). 
[15] B. Utter & R. P. Behringer, Phys. Rev. E 69, 031308 (2004). 
[16] F. Radjai & S. Roux, Phys. Rev. Lett. 92, 174301 (2004). 
[17] J. Schafer, S. Dippel & D. E. Wolf, J. Phys. I(France) 6, 5 (1996). 

[18] D.C. Rapaport, The art of molecular dynamics simulation, Cambridge University Press, Cam- 
bridge UK 2004. 
[19] C. S. Campbell, J. Fluid Mech. 348, 85 (1997). 

[20] R. Arevalo , D. Maza k L. Pugnaloni, Phys. Rev. E 74, 0213031 (2006). 
[21] S. F. Edwards, Granular Matter: An Interdisciplinary Approach, 121, Springer, New York, 
1994. 



17 



[22] J. Geng & R. P. Behringer, Phys. Rev. Lett. 93, 238002 (2004). 

[23] J.-P. Bouchaud & A. Georges, Physics Reports 195, 127 (1990). 

[24] D. H. Zanette, Brazilian Journal of Physics 29, 108 (1999). 

[25] Levy Fligths and Related Topics in Physics. M. F. Shlesinger, G. M. Zaslavsky & U. Frish 

Eds., 269, Springer, Berlin Heidelberg New York, 1995. 

[26] A. Compte, Phys. Rev. E 53, 4191 (1996). 

[27] A. R. Plastino k A. Plastino, Physica A 222, 347 (1995). 

[28] A. Compte & D. Jou, J. Phys. A 29, 4321 (1996). 

[29] L. Borland, Phys. Rev. E 57, 6643 (1998). 

[30] A. Upadhyaya, J-P. Rieu, J. Glazier & Y. Sawada, Physica A 293, 549 (2001). 

[31] K. E. Daniels, C. Beck & E. Bodenschatz, Physica A 193, 208 (2004). 

[32] P. Douglas, S. Bergamini & F. Renzoni, Phys. Rev. Lett. 96, 110601 (2006). 

[33] F. Sattin, J. Phys. A: Math. Gen 36, 1583 (2003). 



18 



